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Abstract. The lowest order nonconforming rectangular element in three dimen- 
sions involves 54 degrees of freedom for the stress and 12 degrees of freedom for 
the displacement. With a modest increase in the number of degrees of freedom (24 
for the stress), we obtain a conforming rectangular element for linear elasticity in 
three dimensions. Moreover, unlike the conforming plane rectangular or simplicial 
elements, this element does not involve any vertex degrees of freedom. Second, we 
remark that further low order elements can be constructed by approximating the 
displacement with rigid body motions. This results in a pair of conforming ele- 
ments with 72 degrees of freedom for the stress and 6 degrees of freedom for the 
displacement. 



1. Introduction 

Let Q C be a contractible polygonal domain occupied by a linearly elastic body 
and let L^(f2,M^) be the space of square integrable vector fields. We denote by S 
the space of 3 x 3 symmetric matrix fields and by H{dw, Q, S) the space of square 
integrable symmetric matrix fields with divergence, taken row-wise, square integrable. 
The stress field a and the displacement u are the unknowns in the elasticity equations 

Aa = e{u), diva = /, inQ. 

Here e{u) = {dui/dxj + duj/dxi)ij=i^__,^3 is the symmetric part of the gradient and 
/ encodes the body forces. The comphance tensor A — A{x) : S — )■ S is given, 
bounded and symmetric positive definite uniformly with respect to x e Jl. The pair 
((7, u) is the unique solution of the Hellinger-Reissner variational formulation: Find 
{a,u) e H{div,n,S) x L'^{n,M.^) such that 



(1.1) 



/ {Aa : T + div T ■ u)dx ^0, re if (div, Q, §), 
Jn 

/ div a ■ V dx — / f-vdx, v & L^{Q,,. 
Jn Jn 



where a : t — CTijTij is the Probenius product of a and r and we assume u — 

on dQ. 

For several decades, the existence of stable, conforming mixed elements for elasticity 
was an open problem. Starting with the seminal work of Arnold and Winther, [8, 9], 
for triangular elements, extended to tetrahedral elements in [1, 4, 18], conforming and 
nonconforming mixed finite elements for elasticity on rectangular meshes have been 
constructed by several authors [3, 11, 27, 28, 20, 21, 15]. We refer to [8, 11, 7, 12] 
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and the references therein for alternative approaches to mixed formulations of the 
elasticity equations. In particular, we note advances on mixed finite elements for 
elasticity where the symmetry of the stress field is enforced weakly using Lagrange 
multipliers, [2, 5, 24, 26, 25, 22, 6, 7, 14, 19, 17, 16, 10]. In general, it requires a great 
deal of ingenuity to construct mixed finite elements for elasticity. 

The element described in this paper is the only one known, in three dimension on 
rectangular meshes, which has strong symmetry for the stress field and is if(div) 
conforming. The tetrahedral analogues have a high number of degrees of freedom 
and involve vertex degrees of freedom. In general, on an arbitrary triangulation, 
vertex degrees of freedom are unavoidable for conforming elements for the stress 
field with strong symmetry [8, 4]. Although it could be expected that on some 
triangulation, elements can be constructed with no vertex degrees of freedon [17], it 
was a surprise to us that this is possible on rectangular meshes in three dimensions 
which are prefered over tetrahedral elements by some engineers. Progress towards 
the discovery of the element was made possible by the observation in [15] that for 
plane elasticity, no vertex degrees of freedom arc necessary for the normal normal 
component of the stress field. Although the resulting element is relatively complicated 
to use, 72 degrees of freedom for the stress or 12 degrees of freedom on average for 
each component of the stress field, we believe the strong symmetry property of the 
element justifies the additional degrees of freedom (the lowest order nonconforming 
element involves 54 degress of freedom for the stress). An example of a preference 
in some situations for elements with strong symmetry for two dimensional elasticity 
is [23] p. 4. We quote "Unfortunately, there exists no simple low order conforming 
element for ..." {H{div,Q,S)). "Possible remedies are to use macro-elements or 
to weaken the symmetry condition, ..." "However, if we weaken the symmetry, an 
additional factor due to the nonconformity of the approach occurs" . We claim that 
for the three dimensional problem, the element described in this paper is the only 
reasonable option available. 

The paper is organized as follows: In the next section we describe the finite element 
spaces on each rectangle after outlining a general framework for the construction of 
mixed finite elements for elasticity. We show how one may obtain shape functions for 
three dimensional rectangular elements for elasticity from the corresponding plane 
elements using tensor products. In Section 3 we prove stability and convergence of 
the mixed finite elements using arguments now standard. We conclude with some 
remarks on low order elements using the space of rigid body motions to approximate 
the displacement and the extension to higher order elements. 



2. Description of the finite elements 

We denote by H^{K,X) the space of functions with domain K C M^, taking values 
in the finite dimensional space X, and with all derivatives of order at most k square 
integrable. For our purposes, X will be either S, ]R^,]R^ or M, and in the latter case, 
we simply write H^{X). The norms in H^{K^ X) and H^{K) are denoted respectively 
by II • IIh* and || ■ ||fc. Let Th denote a conforming partition of Vt into rectangles K 
of size bounded by /i, which is quasi-uniform in the sense that the aspect ratio of the 
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rectangles is bounded by a fixed constant. Whether K is a, plane rectangle or a three 
dimensional brick will be clear from the context. 

We first outline a set of general criteria to satisfy in order to obtain stable mixed finite 
elements for the linear elasticity problem. A conforming mixed finite element approx- 
imation of (1.1) consists in choosing spaces E/j C i?(div, Q,S) and Vh C L^(Q,' 
and determining {ah,Uh) G T^h x Vh such that: 



(2.1) 



{Aah ■ Th + div Th ■ Uh) dx = 0, The T^h, 

/ div ah ■ Vhdx = / / • Vh dx, Vh e 14- 

Jo, Jo, 



We recall that a matrix field r G E/i if and only if rn is continuous across the internal 
faces of Th where n denotes one of the unit vectors normal to the internal face. It is 
now known that for stability, [8, 9, 11, 4], the following two conditions are sufficient: 

• div Eft C Vh 

• There exists a linear operator Uh : H^{fl,Ei) — )■ S/j, such that there is a 
constant c independent of h with ||n/iT||o < c||t||i for all r e H^{fl,E>), and 
such that divn?,r = PhdivT for all r e //^(Q,§) where Ph : L'^{n,R^) Vh 
denotes the projection. 

If we assume that an interpolation operator 11^ : H^{Q, S) — )■ S/^ can be defined which 
reproduces the degrees of freedom of r G T^h, the second condition is satisfied if the 
degrees of freedom for the stress space are chosen such that 

(2.2) / {div UhT - div t) ■ V dx ^ - {UhT - t) : e{v) dx + / {UhT - T)n ■ v ds, 

Jk Jk JdK 

vanishes. Therefore in addition to require normal continuity of the degrees of freedom, 
it is sufficient that the degrees of freedom include j^Tn-v ds and r : e{v), for each 
face / of an element K and each v e Vk- If in addition, the interpolation operator 
reproduces piecewise constant matrix fields, by a standard scaling argument, it is 
bounded and the error analysis also follows from now classical arguments. 

Before describing the new three dimensional conforming elements, we first note that 
the shape functions for the low order three dimensional nonconforming element in 
[21], can be obtained from their two dimensional analogue by a tensor product con- 
struction. 



2.1. Tensor product construction. Wc use the usual notation of Vk{K,X) for 
the space of polynomials on K with values in X of total degree less than k and 
Vki,k2{^j X) for the space of polynomials of degree at most ki in xi and of degree at 
most k2 in X2. Similarly, Vkj^^k2,k3i^^-^) denotes the space of polynomials of degree 
at most ki in xi, of degree at most k2 in X2 and of degree at most ks in x^. We 
write Vk, PkiM ^-iid Vki,k2,k3 respectively when X = M. When the domain K in 
is obvious, Vk{xi,Xj) and Vki^k2{xi,^j) denote respectively the space of polynomials 
of total degree k in the variables Xi and Xj and of degree at most ki in Xi and at 
most k2 in Xj for i not necessarily smaller than j. We also use a similar notation for 
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spaces of polynomials on a face. We will denote by I ^ ^' ^ 7? ^' * J space of 

\' ksM ' ksM/s 

symmetric matrix fields cr — ( '^^^ ^^"^ ) such that uw e Vk^ kii C12 = C2 1 G Vks k^ 

\(T2l (T22J 

and (T22 e Vkr„ke- Similarly the notation [^''^''^^j for Pfci,^ ^ Vk3,k4 will also be used. 

(T^ki.k'zM', '^*:4,fc5,fc6 T^k7,ks,kg \ 

'PkiMM '^kwMiMi 'PkiaMiMr, we mean the space of symmetric matrix fields 
'Pk7,ks,kg '^A;i3,fci4,fci5 'Pkie,kn,kis/ § 

((Jii (7i2 (TlsX 
(721 (722 (723 j such that (7ii e VkiMM^ ^22 £ Vkr,ksM^ ^33 £ 'PkisMtMb^ 
0"31 0"32 (733/ 

(T12 = (721 G Vk^,k5,ke, (7i3 = (731 e VkrMM ^23 = (732 ^ Vk^^MiMs' We will also 

('PkiMM \ 
Vk^MM 1 in place of Vk^MM x 'PkiMM ^ 'PkrMM- Otten at the place of the 

spaces Vk,i,m{K, X), we will write explicit vector spaces as span. 

For two finite dimensional spaces E and F, we denote hj E^F the vector space equal 
to span{e/}, e & E, f & F. A low order plane nonconforming element described in 
uses as shape functions on an element K for the stress 

aeJ:^^(K)=( ^^-^ span{l,xi,X2,a;?,a;|}\ 

"""^^ ' \^span{l,a;i,a;2,a;?,a;|} Pi,i ) 



and for the displacement v e V'^^{K) = ( g^!|^||| j' j ) • postulate that for the 



\^span{l, xi}^ 

three dimensional generalization, the choice of finite dimensional spaces for each com- 

(CTU (7l2 (7l3\ 
(7i2 (722 (723 1 ^ud cach Component of the displacement 
(7l3 (723 CTssJ 

V2 should conform to the following requirements: 



(7ii (7i2 \ ivi \ i^gj^ ^YiQ plane elements in the variables a;i,a;2 

0-12 0-22y \V2j 

when restricted to the planes X3 = 0, 1, 
• I "^^^ '^^^ I and I I belong to the plane elements in the variables xi,xs 
when restricted to the planes X2 — 0, 1, 

(722 (723 I ( ^2 \ 

and belong to the plane elements in the variables X2,X3 

(723 (733/ V^3/ 



when restricted to the planes xi — 0, 1. 
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We note that the above statements can be summarized as follows: for i < j,i,j G 
{1,2,3}, k ^ i'^^^ '^^^ ] and ( ) belong to the plane elements in the vari- 



ables Xi,Xj when restricted to the planes — 0,1, that is an e V\^i{xi,Xj),(Tij e 
span{l, Xj, Xj, xf, and Vi e span{l, Xj}, Vj e span{l, Xi} on the planes Xk — 0, 1. 

A set of shape functions for the three dimensional elasticity problem is obtained by 
taking as finite element spaces for the components of the stress and the displace- 
ment, tensor products of the plane finite element spaces in the variables Xi,Xj with 
span{l,a;fc}. For this note that for i < G {1,2,3}, k ^ 

Vi,i{xi,Xj) ® span{l,Xk} = Pi,i,i, 

span{l, Xi, Xj, x1, x"^} ® span{l, Xk} = span{l, Xi,Xj, x1, x"^, x^, 

2 2 1 

X^Xk, XjXk, X^Xk, XjX^ J 

span{l,a;j} (g) span{l,a;fc}, = span{l, Xi,Xk,XiXk} = Vi^i{xi,Xk). 

We thus obtain the shape functions for the low order three dimensional element 
described in [21]. That is if we denote by S^^f(ii') the finite element for the stress 
and V^^{K) the finite element space for the displacement, a e T,^f{K) if and only 
if an e Pi, 1,1, i = 1,2,3, aij e span{l, Xj, x^-, , x^, XiXk,XjXk,x'lxk,x'^jXk} and 
V e V^^{K) if and only if Vi G Vi,i{xj,Xk) for i < e {1,2, 3}, A; ^ {i,j}. 

Alternatively, if we define 

^2 (-^^i, -^j) span{l, Xi, Xj , Xi iXj, x^, XiX^, XjX^ , Xi x^, XjX^}, 

S^cf(^) = I Vjix^^x^) ^ Vni V2{x2,Xs) 
(P2ixi,Xs) V2{X2,X3) Viu 

2.2. Description of the elements. We first describe the shape functions of a con- 
forming rectangular element recently introduced in [15] . Then, following the strategy 
described in Section 2.1, we give a set of shape functions for the finite element spaces 
in three dimensions. We give a more concise description of the stress finite element 
space before giving a unisolvent set of degrees of freedom. Let us denote by S^P(i^) 
the space of shape functions for the stress for the plane conforming element and by 
S^P(A') the corresponding three dimensional element. Then [15], 

f / 1 2 \/13 2\/2 2 

-,2D ^ ^^20 /" _^ onr>Tn J I 2"^^ X1X2 ^ j 3*^1 X^X2 ^ j X-^^X2 X1X2 



then 



^tr(K) = E-(X) + span^ ( y^2 , ( L ,2 , , . ,2 _ . 3 



/^2 ~,2 



3^1^2/ 

We choose as finite element space for the displacement the same space V^^{K) used in 
the nonconforming case. Note that it can equivalently be written as T'oii x ^loi x ^iio 
and will be denoted by V{K) in the remaining part of this paper. We define T?^[K) 
as the direct sum of Ti^^^{K) and the span of 24 additional shape functions a'^^'^''\i < 
j, i,j = 1,2,3 and A; = 1, . . . , 8 obtained from the additional shape functions in the 
definition of S^P(Ar) in the variables Xi and Xj and the tensor product construction 
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outlined in Section 2.1. For a given matrix, we will refer to the entry on the ith row 
and jth column as the component. 



rl2,(l) 



rl2,(3) 



1^2 


XIX2 




X1X2 




S) 










3-^1 


xfx2 


0^ 


xlx2 


-xixj 















— ^xlxs X1X2X3 0'' 

rl2,(2) _ I X1X2X3 -^X^Xs 

0, 



a 



^xfx^ x\x2Xz 0^ 

12,(4) _ I 2,2^22:3 —Xix\xz 











0> 



(—x\x2 X1X2 0\ I —x\x2X2, XiX^X-i 

xixi -\xlxi , (7^^'(^) = xxxlxz -\xlxz 
0/ \ 0, 

(—'^X\X2 x\x\ 0\ /—f^^'i 3^22^3 x\x\x^ 0^ 

a;?a;| -fxix^ , cr^^'^^^ = x^x^xs -|xix|x3 
oy \ 

Note that the third column of all the above matrices is the zero vector. We conclude 
that for A; = 1,...,8, a^^^fc) ^.qs its (1,3) and (2,3) component 0, but the (1,2) 
component non-zero. 

= I I , a^3,(2) ^ 

(7^3.(5) = 
^13,(7) ^ 



1^2 





XlX3\ 








\,c 


XlX^ 


- 


"2^3/ 


3-^1 





xjxs \ 










x\xz 


- 


-xxxl) 


-xlx3 





X\xX\ 








, 


,xi,x| 





-¥1} 


-|a;5a;3 













•^l'^3 





2 3 
3X1X3 




2 2 
%Xj 1 tX/ 2»*^ 



x\x2x\ — |a;iX2.X3 



Note that the second row or second column of all the matrices cr^^'^'^), k = 1, . . . , 8 
is the zero vector. We conclude that for A; = 1, . . . , 8, a^^'^''^ has its (1, 2) and (2, 3) 
component 0, but the (1,3) component non-zero. 



a 



\ /O 

23,(1) _ I n _la;2 a;2a:3 , a^^'^^^ = -^xix^, X1X2X3 

X2X3 \0 X1X2X3 -\xixl^ 

\ /O 

^23.(3) ^ I -^xl xlx3 , (7'=^'(^' = -\x^xl x^xlx3 

X2X2, —X2X\) \0 X1X2X2, — X1X2X3 
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\ /O 

^23,(5) ^ j Q _3,2^g ^^^2 ^ ^23,(6) ^ q _;5;^a.2^g XiX2xl 
X2X\ ~3^3/ \0 XiX2x\ 

/O \ /O 

(j23,(7) ^ j Q _|a;32,g a;2^2 ^ ^^23,(8) ^ j q -Ix^xlx^ Xixjxj 

\0 ^2^3 — 3*^2*^3/ yO XiX2^3 — 3'^l'^'2'^3y 

Note that the first row or first column of all the matrices cr^^'^'^^ A; = 1, . . . , 8 is the zero 
vector. We conclude that for A; = 1, . . . , 8, cr^^'*^'^) has its (1, 2) and (1, 3) component 
0, but the (1,3) component non-zero. 

We can summarize the above observations by stating that 

(2.3) a'Jf\ A; = 1, . . . , 8 = Ofori < j unless (p, q) = (ij). 

Note that (2.3) can also be seen directly from the tensor product construction. We also 

note that like their two dimensional analogue, the above additional shape functions 
are divergence free and that the dimension of H'^f'^K) is 54 + 24 = 78. This is implicit 
in the tensor product construction and also follows from Lemma 2.1 below. We can 
describe T,lf{K) in a more concise form as follows: r e T,lf{K) if and only if 



(2.4) T e E(X) := \ V221 Vi3i V122 ,divT e 



^311 


^221 


V, 


^221 


7^131 


V 


^212 


"^122 


V 





Clearly Y?^{K') is contained in the latter space. Notice that the dimension of 5*2 

is 102. And for r e ^2, divr G I V\2\ I which has dimension 36. 



P311 


'^221 


V, 


"^221 


Pl31 


V 


^212 


V\22 


V 




Hence the constraint divr G ViK) imposes 36-12=24 constraints. It follows that the 
dimension of S(ir) is at least 78. That its dimension is 78 follows from the unisolvency 
of the degrees of freedom we give in Lemma 2.1 below. Hence, S(i^) = Y^^iK). 

Before we give the degrees of freedom of Y?^{K\ we explain why unlike the con- 
forming plane rectangular or simplicial elements, this element does not involve any 
vertex degrees of freedom. In general, on an arbitrary simplicial triangulation, vertex 
degrees of freedom are necessary [8, 4]. It will be enough to consider the unit cube 
[0, 1]^. We first write exphcitly the condition of normal continuity of 

' cri2 cri3^ 
C12 cr22 cr23 
^VA C23 cr33y 

Clearly on horizontal faces 2:3 = 0, 1, an = ±((Ti3, (T23, o"33)- On vertical faces Xi = 0, 1, 
an — ±(crii, ai2, ai^) and on lateral faces 0:2 = 0, 1 crn = ±(cri2, cr22, 1723). We therefore 
need (733 continuous on faces x^ = 0, 1, a22 continuous on faces a;2 = 0, 1 and an 
continuous on faces Xi — 0, 1. Since for i — 1, 2, 3 the faces Xi — and = 1 do 
not meet, no edge or vertex degrees of freedom are necessary for an. Next, we need 
(T23 continuous on faces x^ = 0, 1, ,T2 = 0, 1, so we need (T23 continuous on the edges 
X2 = 0, X3 = 0; X2 = 1, 2:3 = 0; X2 = 0, X3 = 1; X2 = 1, = 1. And similar conditions 
for (T13 and ai2, that is for i < G {1, 2, 3}, A; ^ we need aij continuous on 
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the faces Xj = 0, 1 and Xj = 0, 1. This requires the continuity of on the four edges 
Xi,Xj G {0, 1}. No vertex degrees of freedom are therefore necessary here. 

Let A' = [0, 1]"^ be the reference cube and let F : K ^ K he an affine mapping onto 
K, F{x) = Bx + b, with b = 62, ^3)^ and 

/hi 
5=0 h2 
\0 hs^ 

To maintain the analogy with the simphcial case, although B is symmetric we will 
write B^. Given a matrix field f : K ^ define r : ^ S by the matrix Piola 
transform 

t{x) = Bt{x)B^. 

It is not difficult to verify that r G S(i^) if and only if f G S(i^). We recall that an 
edge of K is given by Xi,Xj G {0, 1} and we denote by / an integral over 

J{xi=0,l, Xj=0,l} 

such an edge. We now give the degrees of freedom on ^{K). 

Lemma 2.1. An element a G ^{K) is uniquely determined by the following degrees 
of freedom 



(1) / aijV dxf, V G span{l,Xk,} for k ^ i < j i,j = 1,2,3, 

J{x,=0,l, Xj=0,l} 

(24 degrees of freedom), 

(2) / aiivdx,v G 7^i,i(a;;, a;^), i ^ {l,k},i = 1,2,3, (24 degrees of freedom), 

Jxi=0,l 

(3) / aijV dxf, V G span{l,Xk}, k ^ {hj},i < j, I — hj, (24 degrees of 

Jxi=0,l 

freedom), 

(4) / o-ijV dx, V G span{l,Xk}, k ^ {i,j},i < j,i,j = 1,2,3, (6 degrees of 
Jk 

freedom). 

Proof. We assume that all degrees of freedom vanish. Note that on the faces Xi — 0,1, 

o-ij e V2,i{xj,Xk), k ^ i < j,i,j = 1,2,3, (e.g. on the faces X2 = 0,1,0-23 e 

^2,1(3^3, Xi)) and on an edge Xi, xj G {0, 1}, aij G Vi{xk). By the first set of degrees of 
freedom, we have aij = on such an edge. Thus on a face Xi = 0, 1, (jy = Xj{l — Xj)aij, 

Gij G VQ,i{xj,Xk). Since / aijvdxf = 0, i) G span{l,Xjt} by the third set of 

Jxi=0,l 

degrees of freedom, aij = on the faces Xi = 0, 1. Similarly, on the faces xj — 0, 1, 
^ij £ 'P2,i{xi,Xk), k ^ i < = 1,2,3 and by the same argument, &ij = as 

well on the faces xj = 0, 1. 

We therefore have dij = Xj{l — Xj)xi{l — Xi)aij. Since (Tjj G ■P2,2,i(a;i, a;^, x^), we have 
dij G ^0,0,1- 

With the last set of degrees of freedom, we conclude that Uij = 0, i < j, i,j — 1, 2, 3. 



Let us write 

k=l 

a*^'^'^) e M,i < = l,2,3and/c = l,...,8,a- e T.lf{K), 

and recall that 0-12 = 0-13 = 0-23 = 0. We show that q;^^'('=) = q;^^'^*') = a^^'^''^ ^0,k ^ 
1, ... ,8. This is a consequence of the tensor product construction as the following 
argument shows. By (2.3), for i < j 

k=l 

Explicitly 

= ai2 + a^^'^^^TiXa + a^^'^^^xiX2X3 + a^^'^^''xlx2 + a^'^'^^^ x\x2Xz + a^'^'^"'' xyx\ + a^'^'^^'^ x^x^x^ 
+ a^^'^^^x^x^ + Q^'^'^^^ x\x\xj,, 
^13 = cri3 + a^^'^^^xixg + a^^'^^^xixaojg + a^^'^^^XiOTa + a^'^'^^^ x\x2X-i + a^^'^^^ariajg + a^^'^^^aJiaJaXg 

0-23 = 5'23 + oF"^^^X2Xj, + a^^'^^)xiX2a;3 + a^^'^^)a;2X3 + oF'^^'^^ Xxx\x2, + Q;^^'(^)a;2X3 + o?'^'^^^ XxX2x\ 
+ a^^'(^)a;^x^ + a23,(8)^^^2^2^ 

and recall that by definition 

0-12 e Span{l,Xi,X2,Xi,X2,X3,XiX3,X2X3,XiX3,X2X3}, 
0-13 e Span{l,Xi,X3,Xi,X3,X2,Xia:2,X2X3,XiX2,X2X3}, 
0-23 e Span{l,X2,X3,X2,X3,Xi,Xia:2,XiX3,XiX2,XiX3}, 

and the result follows since each of the monomial appearing in the expansion of 
(Tjj, i < j, appears only one time and for each expansion, they are all different. 

We conclude that a e S^^f(ii'). By the second set of degrees of freedom, since 
an e "Pill, an = on each face Xi — 0, 1. Therefore, an — Xi{l — Xijan for some 
polynomial an so (Tjj = for all i = 1, 2, 3. □ 

If we denote by n^,/ = 1,2,3 the normal in the Zth direction, each face is given by 
= 0, 1 for some i G {1,2, 3}, an edge is the intersection of two faces = 0, 1, = 
0, 1 for some i, j. We can then rewrite the degrees of freedom as follows 

(1') arii ■ rij v dx, v e Vi{e), for each edge with normals rii and rij, i ^ j, i,j — 
1,2,3, 

(2') J^an ■ nv dx,v e Vi^i{f), 

(3') J^arii ■ riivdx,, v G span{l,a;fc}, k ^ for each face /, with normal Uj 

and i ^ I and all Z = 1, 2, 3. 
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(0 spanjljXs} span{l,a;2}\ 

span{l,,T3} span{l,Xi} , where 

span{l,a;2} span{l,a;i} 
the zero vector space is simply denoted 0. We can thus write the fourth set of degrees 
of freedom as 

(4') f^a:e{v), veV{K). 

Given the discussion before the above lemma, the degrees of freedom in Lemma 2.1 
imply that for a e "^{K), an is uniquely determined by the degrees of freedom on 
the face with normal n. 

On the rectangular partition Th, we define S/j as the space of matrix fields which 
belong piecewise to Ti{K) subject to the continuity conditions that rn for r G S/j is 
continuous across faces. We denote by Vh the space of vector fields which belongs to 
Vk for each rectangle K e Th- 



3. Stability and Error analysis 

Our finite element spaces satisfy E/^ C H{div, S) and Vh C L'^{fl, M) with div T,h C Vh. 
We next define an interpolation operator Uh : H^{fl, §) — )■ T,h with the required prop- 
erties. Because of the edge degrees of freedom, the canonical interpolation operator 
for E/j is not well-defined on H^(fl,E>). As in [3], we first define an interpolation 
operator : H\k,E,) E{k) by 

J{Il^^{a))ni ■ rij i) dx = e 7^i(e) for e with normals andn^, i ^ j, i,j = 1, 2, 3, 

J^{a-U%{a))n-nvdx^O "iv eVi,i{f), 

a — n^(a"))ni • nidx — Q for each face /, with normal rii and i ^ I and alH = 1, 2, 3, 

[ {a- U%{a)) ■.€{{;) ^0 \/v eV{k). 
Jk 

By Lemma 2.1, B^ is well defined. Next define U% : H\K,E,) ^K) by 

U.%t{x) = BU^^f{x)B^, 
for each rectangle K oiTh and define B° : H^{Q, S) — > E/j by 

(B°r)U = B^T. 
A standard scaling argument gives 
(3.1) ||B°r-r||o<c/i||r||i, 

where c does not depend on h since B° reproduces piecewise constant matrix fields. 
It follows that B° is bounded on H^{Q,S). 

The proof of the commutativity property divB^r = P/jdivr for all r G H^{Q,Ei) is 
similar to the one given in [3], pp. 1423-1424. Next, let Rh be a Clement interpolation 
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operator [13] which maps L'^{Q, S) into the subspace of E/j of continuous matrix fields 
whose components are piecewise in -Pi, 1,1. We have 

WRhT - T\\j < ch'^-mrWrn, < J < 1, J < m < 2, 
with c independent of h. The interpolation operator is defined by 

(3.2) Uh^Ul{I-Rh) + Rh, 

and as in [3], one establishes its boundedness and the commutativity property. We 
conclude that the pair (E/i, Vh) is stable. 

For the error analysis, we recall that the projection operator satisfies the error 
estimate 

(3.3) \\PhV -v\\o< ch'^WvWm, ior aWv e H"\n), < m < 1, 

since \4 contains piecewise constant vector fields. We then have the following theorem 

Theorem 3.1. Let {a,u) and {ah,Uh) be the unique solutions of (1.1) and (2.1) 
respectively. Then 

II divo" — div(Th||o < c/i™|| divcrU^, < m < 1, 
\\o- - (ThWo < ch\\a\\i, 
\\u — Uh\\o < ch\\u\\2- 

Proof. We only outhne the proof since it is essentially the same as the one in [8]. 
Define ||.|U by ||t||^ = /n^''' ■ "^dx. Then one shows that 

Ik - cthWa < \\cr - ^hCr\\A- 
Since the norm ||.||^ is equivalent to the norm, we obtain using (3.1), 

\\0- - CTfello < c/l||o-||i. 

One also shows that divert = P/j diver. The commutativity property and (3.3) yields 
the error estimate on diver. For the error estimate for the displacement, again from 
[8] , we have 

\\PhU - Uh\\o <\\cr- (Th\\o, 

so 

\\u - Uh\\o < \\u - Phu\\o + \\PhU - Uh\\o < ch\\u\\i + ch\\a\\i < ch\\u\\2, 
since Aa — €{u) . □ 

4. Concluding remarks 

In three dimensions, the space of rigid body motions T is the space of linear polynomial 
functions of the form x i-^ a + b x x for some a,b & M:'. It is contained in V{K) for 
each clement K. Hence, as in [8, 9, 3, 4, 11], we can take V{K) — T and the stress 
space as 



(4.1) TeJ:{K):={ \V221 ^131 ^122 ,divreT 



^311 


^221 


V, 


^221 


^131 


V 


^212 


^122 


V 
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Since for v G T,e(f) = 0, the degrees of freedom are the same as the ones in Lemma 
2.1 with the last set of degrees of freedom removed. The error analysis is the same. 
We note that the elements of [20, 21, 15] can be simplified as well this way. 

As noted for example in [12], high order mixed finite elements for elasticity usually 
require a high number of degrees of freedom and typically require to find the dimension 
and a base for the space 

{r e E(ir), diver = 0, crn = 0, ondQ}. 

A possible approach is to follow the proof in the simplicial case [4] , but this goes well 
beyond the scope of these remarks. Thus for that reason as well, we did not pursue 
here the construction of a discrete elasticity sequence. 
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